<!DOCTYPE html>

<html xmlns="http://www.w3.org/1999/xhtml">

<head>

<meta charset="utf-8">
<meta http-equiv="Content-Type" content="text/html; charset=utf-8" />
<meta name="generator" content="pandoc" />
<meta name="viewport" content="width=device-width, initial-scale=1">

<style type="text/css">
@font-face {
font-family: octicons-link;
src: url(data:font/woff;charset=utf-8;base64,d09GRgABAAAAAAZwABAAAAAACFQAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAABEU0lHAAAGaAAAAAgAAAAIAAAAAUdTVUIAAAZcAAAACgAAAAoAAQAAT1MvMgAAAyQAAABJAAAAYFYEU3RjbWFwAAADcAAAAEUAAACAAJThvmN2dCAAAATkAAAABAAAAAQAAAAAZnBnbQAAA7gAAACyAAABCUM+8IhnYXNwAAAGTAAAABAAAAAQABoAI2dseWYAAAFsAAABPAAAAZwcEq9taGVhZAAAAsgAAAA0AAAANgh4a91oaGVhAAADCAAAABoAAAAkCA8DRGhtdHgAAAL8AAAADAAAAAwGAACfbG9jYQAAAsAAAAAIAAAACABiATBtYXhwAAACqAAAABgAAAAgAA8ASm5hbWUAAAToAAABQgAAAlXu73sOcG9zdAAABiwAAAAeAAAAME3QpOBwcmVwAAAEbAAAAHYAAAB/aFGpk3jaTY6xa8JAGMW/O62BDi0tJLYQincXEypYIiGJjSgHniQ6umTsUEyLm5BV6NDBP8Tpts6F0v+k/0an2i+itHDw3v2+9+DBKTzsJNnWJNTgHEy4BgG3EMI9DCEDOGEXzDADU5hBKMIgNPZqoD3SilVaXZCER3/I7AtxEJLtzzuZfI+VVkprxTlXShWKb3TBecG11rwoNlmmn1P2WYcJczl32etSpKnziC7lQyWe1smVPy/Lt7Kc+0vWY/gAgIIEqAN9we0pwKXreiMasxvabDQMM4riO+qxM2ogwDGOZTXxwxDiycQIcoYFBLj5K3EIaSctAq2kTYiw+ymhce7vwM9jSqO8JyVd5RH9gyTt2+J/yUmYlIR0s04n6+7Vm1ozezUeLEaUjhaDSuXHwVRgvLJn1tQ7xiuVv/ocTRF42mNgZGBgYGbwZOBiAAFGJBIMAAizAFoAAABiAGIAznjaY2BkYGAA4in8zwXi+W2+MjCzMIDApSwvXzC97Z4Ig8N/BxYGZgcgl52BCSQKAA3jCV8CAABfAAAAAAQAAEB42mNgZGBg4f3vACQZQABIMjKgAmYAKEgBXgAAeNpjYGY6wTiBgZWBg2kmUxoDA4MPhGZMYzBi1AHygVLYQUCaawqDA4PChxhmh/8ODDEsvAwHgMKMIDnGL0x7gJQCAwMAJd4MFwAAAHjaY2BgYGaA4DAGRgYQkAHyGMF8NgYrIM3JIAGVYYDT+AEjAwuDFpBmA9KMDEwMCh9i/v8H8sH0/4dQc1iAmAkALaUKLgAAAHjaTY9LDsIgEIbtgqHUPpDi3gPoBVyRTmTddOmqTXThEXqrob2gQ1FjwpDvfwCBdmdXC5AVKFu3e5MfNFJ29KTQT48Ob9/lqYwOGZxeUelN2U2R6+cArgtCJpauW7UQBqnFkUsjAY/kOU1cP+DAgvxwn1chZDwUbd6CFimGXwzwF6tPbFIcjEl+vvmM/byA48e6tWrKArm4ZJlCbdsrxksL1AwWn/yBSJKpYbq8AXaaTb8AAHja28jAwOC00ZrBeQNDQOWO//sdBBgYGRiYWYAEELEwMTE4uzo5Zzo5b2BxdnFOcALxNjA6b2ByTswC8jYwg0VlNuoCTWAMqNzMzsoK1rEhNqByEyerg5PMJlYuVueETKcd/89uBpnpvIEVomeHLoMsAAe1Id4AAAAAAAB42oWQT07CQBTGv0JBhagk7HQzKxca2sJCE1hDt4QF+9JOS0nbaaYDCQfwCJ7Au3AHj+LO13FMmm6cl7785vven0kBjHCBhfpYuNa5Ph1c0e2Xu3jEvWG7UdPDLZ4N92nOm+EBXuAbHmIMSRMs+4aUEd4Nd3CHD8NdvOLTsA2GL8M9PODbcL+hD7C1xoaHeLJSEao0FEW14ckxC+TU8TxvsY6X0eLPmRhry2WVioLpkrbp84LLQPGI7c6sOiUzpWIWS5GzlSgUzzLBSikOPFTOXqly7rqx0Z1Q5BAIoZBSFihQYQOOBEdkCOgXTOHA07HAGjGWiIjaPZNW13/+lm6S9FT7rLHFJ6fQbkATOG1j2OFMucKJJsxIVfQORl+9Jyda6Sl1dUYhSCm1dyClfoeDve4qMYdLEbfqHf3O/AdDumsjAAB42mNgYoAAZQYjBmyAGYQZmdhL8zLdDEydARfoAqIAAAABAAMABwAKABMAB///AA8AAQAAAAAAAAAAAAAAAAABAAAAAA==) format('woff');
}
body {
-webkit-text-size-adjust: 100%;
text-size-adjust: 100%;
color: #333;
font-family: "Helvetica Neue", Helvetica, "Segoe UI", Arial, freesans, sans-serif, "Apple Color Emoji", "Segoe UI Emoji", "Segoe UI Symbol";
font-size: 16px;
line-height: 1.6;
word-wrap: break-word;
}
a {
background-color: transparent;
}
a:active,
a:hover {
outline: 0;
}
strong {
font-weight: bold;
}
h1 {
font-size: 2em;
margin: 0.67em 0;
}
img {
border: 0;
}
hr {
box-sizing: content-box;
height: 0;
}
pre {
overflow: auto;
}
code,
kbd,
pre {
font-family: monospace, monospace;
font-size: 1em;
}
input {
color: inherit;
font: inherit;
margin: 0;
}
html input[disabled] {
cursor: default;
}
input {
line-height: normal;
}
input[type="checkbox"] {
box-sizing: border-box;
padding: 0;
}
table {
border-collapse: collapse;
border-spacing: 0;
}
td,
th {
padding: 0;
}
* {
box-sizing: border-box;
}
input {
font: 13px / 1.4 Helvetica, arial, nimbussansl, liberationsans, freesans, clean, sans-serif, "Apple Color Emoji", "Segoe UI Emoji", "Segoe UI Symbol";
}
a {
color: #4078c0;
text-decoration: none;
}
a:hover,
a:active {
text-decoration: underline;
}
hr {
height: 0;
margin: 15px 0;
overflow: hidden;
background: transparent;
border: 0;
border-bottom: 1px solid #ddd;
}
hr:before {
display: table;
content: "";
}
hr:after {
display: table;
clear: both;
content: "";
}
h1,
h2,
h3,
h4,
h5,
h6 {
margin-top: 15px;
margin-bottom: 15px;
line-height: 1.1;
}
h1 {
font-size: 30px;
}
h2 {
font-size: 21px;
}
h3 {
font-size: 16px;
}
h4 {
font-size: 14px;
}
h5 {
font-size: 12px;
}
h6 {
font-size: 11px;
}
blockquote {
margin: 0;
}
ul,
ol {
padding: 0;
margin-top: 0;
margin-bottom: 0;
}
ol ol,
ul ol {
list-style-type: lower-roman;
}
ul ul ol,
ul ol ol,
ol ul ol,
ol ol ol {
list-style-type: lower-alpha;
}
dd {
margin-left: 0;
}
code {
font-family: Consolas, "Liberation Mono", Menlo, Courier, monospace;
font-size: 12px;
}
pre {
margin-top: 0;
margin-bottom: 0;
font: 12px Consolas, "Liberation Mono", Menlo, Courier, monospace;
}
.select::-ms-expand {
opacity: 0;
}
.octicon {
font: normal normal normal 16px/1 octicons-link;
display: inline-block;
text-decoration: none;
text-rendering: auto;
-webkit-font-smoothing: antialiased;
-moz-osx-font-smoothing: grayscale;
-webkit-user-select: none;
-moz-user-select: none;
-ms-user-select: none;
user-select: none;
}
.octicon-link:before {
content: '\f05c';
}
.markdown-body:before {
display: table;
content: "";
}
.markdown-body:after {
display: table;
clear: both;
content: "";
}
.markdown-body>*:first-child {
margin-top: 0 !important;
}
.markdown-body>*:last-child {
margin-bottom: 0 !important;
}
a:not([href]) {
color: inherit;
text-decoration: none;
}
.anchor {
display: inline-block;
padding-right: 2px;
margin-left: -18px;
}
.anchor:focus {
outline: none;
}
h1,
h2,
h3,
h4,
h5,
h6 {
margin-top: 1em;
margin-bottom: 16px;
font-weight: bold;
line-height: 1.4;
}
h1 .octicon-link,
h2 .octicon-link,
h3 .octicon-link,
h4 .octicon-link,
h5 .octicon-link,
h6 .octicon-link {
color: #000;
vertical-align: middle;
visibility: hidden;
}
h1:hover .anchor,
h2:hover .anchor,
h3:hover .anchor,
h4:hover .anchor,
h5:hover .anchor,
h6:hover .anchor {
text-decoration: none;
}
h1:hover .anchor .octicon-link,
h2:hover .anchor .octicon-link,
h3:hover .anchor .octicon-link,
h4:hover .anchor .octicon-link,
h5:hover .anchor .octicon-link,
h6:hover .anchor .octicon-link {
visibility: visible;
}
h1 {
padding-bottom: 0.3em;
font-size: 2.25em;
line-height: 1.2;
border-bottom: 1px solid #eee;
}
h1 .anchor {
line-height: 1;
}
h2 {
padding-bottom: 0.3em;
font-size: 1.75em;
line-height: 1.225;
border-bottom: 1px solid #eee;
}
h2 .anchor {
line-height: 1;
}
h3 {
font-size: 1.5em;
line-height: 1.43;
}
h3 .anchor {
line-height: 1.2;
}
h4 {
font-size: 1.25em;
}
h4 .anchor {
line-height: 1.2;
}
h5 {
font-size: 1em;
}
h5 .anchor {
line-height: 1.1;
}
h6 {
font-size: 1em;
color: #777;
}
h6 .anchor {
line-height: 1.1;
}
p,
blockquote,
ul,
ol,
dl,
table,
pre {
margin-top: 0;
margin-bottom: 16px;
}
hr {
height: 4px;
padding: 0;
margin: 16px 0;
background-color: #e7e7e7;
border: 0 none;
}
ul,
ol {
padding-left: 2em;
}
ul ul,
ul ol,
ol ol,
ol ul {
margin-top: 0;
margin-bottom: 0;
}
li>p {
margin-top: 16px;
}
dl {
padding: 0;
}
dl dt {
padding: 0;
margin-top: 16px;
font-size: 1em;
font-style: italic;
font-weight: bold;
}
dl dd {
padding: 0 16px;
margin-bottom: 16px;
}
blockquote {
padding: 0 15px;
color: #777;
border-left: 4px solid #ddd;
}
blockquote>:first-child {
margin-top: 0;
}
blockquote>:last-child {
margin-bottom: 0;
}
table {
display: block;
width: 100%;
overflow: auto;
word-break: normal;
word-break: keep-all;
}
table th {
font-weight: bold;
}
table th,
table td {
padding: 6px 13px;
border: 1px solid #ddd;
}
table tr {
background-color: #fff;
border-top: 1px solid #ccc;
}
table tr:nth-child(2n) {
background-color: #f8f8f8;
}
img {
max-width: 100%;
box-sizing: content-box;
background-color: #fff;
}
code {
padding: 0;
padding-top: 0.2em;
padding-bottom: 0.2em;
margin: 0;
font-size: 85%;
background-color: rgba(0,0,0,0.04);
border-radius: 3px;
}
code:before,
code:after {
letter-spacing: -0.2em;
content: "\00a0";
}
pre>code {
padding: 0;
margin: 0;
font-size: 100%;
word-break: normal;
white-space: pre;
background: transparent;
border: 0;
}
.highlight {
margin-bottom: 16px;
}
.highlight pre,
pre {
padding: 16px;
overflow: auto;
font-size: 85%;
line-height: 1.45;
background-color: #f7f7f7;
border-radius: 3px;
}
.highlight pre {
margin-bottom: 0;
word-break: normal;
}
pre {
word-wrap: normal;
}
pre code {
display: inline;
max-width: initial;
padding: 0;
margin: 0;
overflow: initial;
line-height: inherit;
word-wrap: normal;
background-color: transparent;
border: 0;
}
pre code:before,
pre code:after {
content: normal;
}
kbd {
display: inline-block;
padding: 3px 5px;
font-size: 11px;
line-height: 10px;
color: #555;
vertical-align: middle;
background-color: #fcfcfc;
border: solid 1px #ccc;
border-bottom-color: #bbb;
border-radius: 3px;
box-shadow: inset 0 -1px 0 #bbb;
}
.pl-c {
color: #969896;
}
.pl-c1,
.pl-s .pl-v {
color: #0086b3;
}
.pl-e,
.pl-en {
color: #795da3;
}
.pl-s .pl-s1,
.pl-smi {
color: #333;
}
.pl-ent {
color: #63a35c;
}
.pl-k {
color: #a71d5d;
}
.pl-pds,
.pl-s,
.pl-s .pl-pse .pl-s1,
.pl-sr,
.pl-sr .pl-cce,
.pl-sr .pl-sra,
.pl-sr .pl-sre {
color: #183691;
}
.pl-v {
color: #ed6a43;
}
.pl-id {
color: #b52a1d;
}
.pl-ii {
background-color: #b52a1d;
color: #f8f8f8;
}
.pl-sr .pl-cce {
color: #63a35c;
font-weight: bold;
}
.pl-ml {
color: #693a17;
}
.pl-mh,
.pl-mh .pl-en,
.pl-ms {
color: #1d3e81;
font-weight: bold;
}
.pl-mq {
color: #008080;
}
.pl-mi {
color: #333;
font-style: italic;
}
.pl-mb {
color: #333;
font-weight: bold;
}
.pl-md {
background-color: #ffecec;
color: #bd2c00;
}
.pl-mi1 {
background-color: #eaffea;
color: #55a532;
}
.pl-mdr {
color: #795da3;
font-weight: bold;
}
.pl-mo {
color: #1d3e81;
}
kbd {
display: inline-block;
padding: 3px 5px;
font: 11px Consolas, "Liberation Mono", Menlo, Courier, monospace;
line-height: 10px;
color: #555;
vertical-align: middle;
background-color: #fcfcfc;
border: solid 1px #ccc;
border-bottom-color: #bbb;
border-radius: 3px;
box-shadow: inset 0 -1px 0 #bbb;
}
.task-list-item {
list-style-type: none;
}
.task-list-item+.task-list-item {
margin-top: 3px;
}
.task-list-item input {
margin: 0 0.35em 0.25em -1.6em;
vertical-align: middle;
}
:checked+.radio-label {
z-index: 1;
position: relative;
border-color: #4078c0;
}
.sourceLine {
display: inline-block;
}
code .kw { color: #000000; }
code .dt { color: #ed6a43; }
code .dv { color: #009999; }
code .bn { color: #009999; }
code .fl { color: #009999; }
code .ch { color: #009999; }
code .st { color: #183691; }
code .co { color: #969896; }
code .ot { color: #0086b3; }
code .al { color: #a61717; }
code .fu { color: #63a35c; }
code .er { color: #a61717; background-color: #e3d2d2; }
code .wa { color: #000000; }
code .cn { color: #008080; }
code .sc { color: #008080; }
code .vs { color: #183691; }
code .ss { color: #183691; }
code .im { color: #000000; }
code .va {color: #008080; }
code .cf { color: #000000; }
code .op { color: #000000; }
code .bu { color: #000000; }
code .ex { color: #000000; }
code .pp { color: #999999; }
code .at { color: #008080; }
code .do { color: #969896; }
code .an { color: #008080; }
code .cv { color: #008080; }
code .in { color: #008080; }
</style>
<style>
body {
  box-sizing: border-box;
  min-width: 200px;
  max-width: 980px;
  margin: 0 auto;
  padding: 45px;
  padding-top: 0px;
}
</style>

</head>

<body>

<!-- README.md is generated from README.Rmd. Please edit that file -->

<h1 id="rhosa--higher-order-spectral-analysis-in-r">rhosa – Higher-Order Spectral Analysis in R</h1>
<!-- badges: start -->

<p><a href="https://travis-ci.com/tabe/rhosa"><img src="" alt="Build Status" /></a> <a href="https://CRAN.R-project.org/package=rhosa"><img src="data:image/svg+xml; charset=utf-8;base64,PHN2ZyB4bWxucz0iaHR0cDovL3d3dy53My5vcmcvMjAwMC9zdmciIHdpZHRoPSI4NSIgaGVpZ2h0PSIyMCI+CiAgPGxpbmVhckdyYWRpZW50IGlkPSJiIiB4Mj0iMCIgeTI9IjEwMCUiPgogICAgPHN0b3Agb2Zmc2V0PSIwIiBzdG9wLWNvbG9yPSIjYmJiIiBzdG9wLW9wYWNpdHk9Ii4xIi8+CiAgICA8c3RvcCBvZmZzZXQ9IjEiIHN0b3Atb3BhY2l0eT0iLjEiLz4KICA8L2xpbmVhckdyYWRpZW50PgogIDxtYXNrIGlkPSJhIj4KICAgIDxyZWN0IHdpZHRoPSI4NSIgaGVpZ2h0PSIyMCIgcng9IjMiIGZpbGw9IiNmZmYiLz4KICA8L21hc2s+CiAgPGcgbWFzaz0idXJsKCNhKSI+CiAgICA8cGF0aCBmaWxsPSIjNTU1IiBkPSJNMCAwaDQzdjIwSDB6Ii8+CiAgICA8cGF0aCBmaWxsPSIjNGMxIiBkPSJNNDMgMGg2M3YyMEg0M3oiLz4KICAgIDxwYXRoIGZpbGw9InVybCgjYikiIGQ9Ik0wIDBoODV2MjBIMHoiLz4KICA8L2c+CiAgPGcgZmlsbD0iI2ZmZiIgdGV4dC1hbmNob3I9Im1pZGRsZSIKICAgICBmb250LWZhbWlseT0iRGVqYVZ1IFNhbnMsVmVyZGFuYSxHZW5ldmEsc2Fucy1zZXJpZiIgZm9udC1zaXplPSIxMSI+CiAgICA8dGV4dCB4PSIyMS41IiB5PSIxNSIgZmlsbD0iIzAxMDEwMSIgZmlsbC1vcGFjaXR5PSIuMyI+CiAgICAgIENSQU4KICAgIDwvdGV4dD4KICAgIDx0ZXh0IHg9IjIxLjUiIHk9IjE0Ij4KICAgICAgQ1JBTgogICAgPC90ZXh0PgogICAgPHRleHQgeD0iNjMiIHk9IjE1IiBmaWxsPSIjMDEwMTAxIiBmaWxsLW9wYWNpdHk9Ii4zIj4KICAgICAgMC4xLjAKICAgIDwvdGV4dD4KICAgIDx0ZXh0IHg9IjYzIiB5PSIxNCI+CiAgICAgIDAuMS4wCiAgICA8L3RleHQ+CiAgPC9nPgo8L3N2Zz4=" alt="CRAN status" /></a></p>
<!-- badges: end -->

<p>This package aims to provide functions to analyze and estimate higher-order spectra or polyspectra of multivariate time series, such as bispectrum and bicoherence [1].</p>
<h2 id="installation">Installation</h2>
<p>You can install the released version of rhosa from <a href="https://CRAN.R-project.org">CRAN</a> with:</p>
<div class="sourceCode" id="cb1"><pre class="sourceCode r"><code class="sourceCode r"><span id="cb1-1"><a href="#cb1-1"></a><span class="kw">install.packages</span>(<span class="st">&quot;rhosa&quot;</span>)</span></code></pre></div>
<p>Alternatively, the development version from <a href="https://github.com/">GitHub</a> with <a href="https://cran.r-project.org/package=remotes">remotes</a>:</p>
<div class="sourceCode" id="cb2"><pre class="sourceCode r"><code class="sourceCode r"><span id="cb2-1"><a href="#cb2-1"></a><span class="co"># install.packages(&quot;remotes&quot;)</span></span>
<span id="cb2-2"><a href="#cb2-2"></a>remotes<span class="op">::</span><span class="kw">install_github</span>(<span class="st">&quot;tabe/rhosa&quot;</span>)</span></code></pre></div>
<h2 id="example">Example</h2>
<p>This is a simple example, based on the outline at Figure 1 of [2], which demonstrates how to use rhosa’s functions to find an obscure relationship between two frequencies in some time series imitated by a generative model.</p>
<p>With four cosinusoidal waves having arbitrarily different phases (<code>omega_a</code>, <code>omega_b</code>, <code>omega_c</code>, and <code>omega_d</code>), but sharing a couple of frequencies (<code>f_1</code> and <code>f_2</code>), we define function <code>D(t)</code> to simulate a pair of time series: <code>v</code> and <code>w</code>. We make them noisy by adding an independent random variate that follows the standard normal distribution.</p>
<div class="sourceCode" id="cb3"><pre class="sourceCode r"><code class="sourceCode r"><span id="cb3-1"><a href="#cb3-1"></a><span class="kw">set.seed</span>(<span class="dv">1</span>)</span>
<span id="cb3-2"><a href="#cb3-2"></a>f_<span class="dv">1</span> &lt;-<span class="st"> </span><span class="fl">0.35</span></span>
<span id="cb3-3"><a href="#cb3-3"></a>f_<span class="dv">2</span> &lt;-<span class="st"> </span><span class="fl">0.2</span></span>
<span id="cb3-4"><a href="#cb3-4"></a>D &lt;-<span class="st"> </span><span class="cf">function</span>(t) {</span>
<span id="cb3-5"><a href="#cb3-5"></a>    omega_a &lt;-<span class="st"> </span><span class="kw">runif</span>(<span class="dv">1</span>, <span class="dt">min =</span> <span class="dv">0</span>, <span class="dt">max =</span> <span class="dv">2</span> <span class="op">*</span><span class="st"> </span>pi)</span>
<span id="cb3-6"><a href="#cb3-6"></a>    omega_b &lt;-<span class="st"> </span><span class="kw">runif</span>(<span class="dv">1</span>, <span class="dt">min =</span> <span class="dv">0</span>, <span class="dt">max =</span> <span class="dv">2</span> <span class="op">*</span><span class="st"> </span>pi)</span>
<span id="cb3-7"><a href="#cb3-7"></a>    omega_c &lt;-<span class="st"> </span><span class="kw">runif</span>(<span class="dv">1</span>, <span class="dt">min =</span> <span class="dv">0</span>, <span class="dt">max =</span> <span class="dv">2</span> <span class="op">*</span><span class="st"> </span>pi)</span>
<span id="cb3-8"><a href="#cb3-8"></a>    omega_d &lt;-<span class="st"> </span><span class="kw">runif</span>(<span class="dv">1</span>, <span class="dt">min =</span> <span class="dv">0</span>, <span class="dt">max =</span> <span class="dv">2</span> <span class="op">*</span><span class="st"> </span>pi)</span>
<span id="cb3-9"><a href="#cb3-9"></a>    wave_a &lt;-<span class="st"> </span><span class="cf">function</span>(t) <span class="kw">cos</span>(<span class="dv">2</span> <span class="op">*</span><span class="st"> </span>pi <span class="op">*</span><span class="st"> </span>f_<span class="dv">1</span> <span class="op">*</span><span class="st"> </span>t <span class="op">+</span><span class="st"> </span>omega_a)</span>
<span id="cb3-10"><a href="#cb3-10"></a>    wave_b &lt;-<span class="st"> </span><span class="cf">function</span>(t) <span class="kw">cos</span>(<span class="dv">2</span> <span class="op">*</span><span class="st"> </span>pi <span class="op">*</span><span class="st"> </span>f_<span class="dv">2</span> <span class="op">*</span><span class="st"> </span>t <span class="op">+</span><span class="st"> </span>omega_b)</span>
<span id="cb3-11"><a href="#cb3-11"></a>    wave_c &lt;-<span class="st"> </span><span class="cf">function</span>(t) <span class="kw">cos</span>(<span class="dv">2</span> <span class="op">*</span><span class="st"> </span>pi <span class="op">*</span><span class="st"> </span>f_<span class="dv">1</span> <span class="op">*</span><span class="st"> </span>t <span class="op">+</span><span class="st"> </span>omega_c)</span>
<span id="cb3-12"><a href="#cb3-12"></a>    wave_d &lt;-<span class="st"> </span><span class="cf">function</span>(t) <span class="kw">cos</span>(<span class="dv">2</span> <span class="op">*</span><span class="st"> </span>pi <span class="op">*</span><span class="st"> </span>f_<span class="dv">2</span> <span class="op">*</span><span class="st"> </span>t <span class="op">+</span><span class="st"> </span>omega_d)</span>
<span id="cb3-13"><a href="#cb3-13"></a>    curve_v &lt;-<span class="st"> </span><span class="cf">function</span>(t) <span class="kw">wave_a</span>(t) <span class="op">+</span><span class="st"> </span><span class="kw">wave_b</span>(t) <span class="op">+</span><span class="st"> </span><span class="kw">wave_a</span>(t) <span class="op">*</span><span class="st"> </span><span class="kw">wave_b</span>(t)</span>
<span id="cb3-14"><a href="#cb3-14"></a>    curve_w &lt;-<span class="st"> </span><span class="cf">function</span>(t) <span class="kw">wave_c</span>(t) <span class="op">+</span><span class="st"> </span><span class="kw">wave_d</span>(t) <span class="op">+</span><span class="st"> </span><span class="kw">wave_c</span>(t) <span class="op">*</span><span class="st"> </span><span class="kw">wave_b</span>(t)</span>
<span id="cb3-15"><a href="#cb3-15"></a>    <span class="kw">data.frame</span>(<span class="dt">v =</span> <span class="kw">curve_v</span>(t) <span class="op">+</span><span class="st"> </span><span class="kw">rnorm</span>(<span class="kw">length</span>(t)),</span>
<span id="cb3-16"><a href="#cb3-16"></a>               <span class="dt">w =</span> <span class="kw">curve_w</span>(t) <span class="op">+</span><span class="st"> </span><span class="kw">rnorm</span>(<span class="kw">length</span>(t)))</span>
<span id="cb3-17"><a href="#cb3-17"></a>}</span></code></pre></div>
<p>Both <code>v</code> and <code>w</code> are oscillatory in principle:</p>
<div class="sourceCode" id="cb4"><pre class="sourceCode r"><code class="sourceCode r"><span id="cb4-1"><a href="#cb4-1"></a>data &lt;-<span class="st"> </span><span class="kw">D</span>(<span class="kw">seq_len</span>(<span class="dv">2048</span>))</span>
<span id="cb4-2"><a href="#cb4-2"></a><span class="kw">with</span>(data, {</span>
<span id="cb4-3"><a href="#cb4-3"></a>    <span class="kw">plot</span>(<span class="kw">seq_len</span>(<span class="dv">100</span>), <span class="kw">head</span>(v, <span class="dv">100</span>), <span class="dt">type =</span> <span class="st">&quot;l&quot;</span>, <span class="dt">col =</span> <span class="st">&quot;green&quot;</span>, <span class="dt">ylim =</span> <span class="kw">c</span>(<span class="op">-</span><span class="dv">3</span>, <span class="dv">3</span>), <span class="dt">xlab =</span> <span class="st">&quot;t&quot;</span>, <span class="dt">ylab =</span> <span class="st">&quot;value&quot;</span>)</span>
<span id="cb4-4"><a href="#cb4-4"></a>    <span class="kw">lines</span>(<span class="kw">seq_len</span>(<span class="dv">100</span>), <span class="kw">head</span>(w, <span class="dv">100</span>), <span class="dt">col =</span> <span class="st">&quot;orange&quot;</span>)</span>
<span id="cb4-5"><a href="#cb4-5"></a>})</span></code></pre></div>
<div class="figure">

<img src="" alt="v and w." width="50%" />

<p class="caption">

<p>v and w.</p>
</p>

</div>

<p>It is noteworthy that the power spectrum densities of <code>v</code> and <code>w</code> are basically identical as shown in their spectral density estimation:</p>
<div class="sourceCode" id="cb5"><pre class="sourceCode r"><code class="sourceCode r"><span id="cb5-1"><a href="#cb5-1"></a><span class="kw">with</span>(data, {</span>
<span id="cb5-2"><a href="#cb5-2"></a>    <span class="kw">spectrum</span>(v, <span class="dt">main =</span> <span class="st">&quot;v&quot;</span>, <span class="dt">col =</span> <span class="st">&quot;green&quot;</span>)</span>
<span id="cb5-3"><a href="#cb5-3"></a>    <span class="kw">spectrum</span>(w, <span class="dt">main =</span> <span class="st">&quot;w&quot;</span>, <span class="dt">col =</span> <span class="st">&quot;orange&quot;</span>)</span>
<span id="cb5-4"><a href="#cb5-4"></a>})</span></code></pre></div>
<div class="figure">

<p><img src="" alt="Spectral density estimation via periodograms." width="50%" /><img src="" alt="Spectral density estimation via periodograms." width="50%" /></p>
<p class="caption">

<p>Spectral density estimation via periodograms.</p>
</p>

</div>

<p>On the other hand, their bispectra are different. More specifically, we are going to see that their bicoherence at some pairs of frequencies are different. rhosa’s <code>bicoherence</code> function allows us to estimate the magnitude-squared bicoherence from samples.</p>
<div class="sourceCode" id="cb6"><pre class="sourceCode r"><code class="sourceCode r"><span id="cb6-1"><a href="#cb6-1"></a>x &lt;-<span class="st"> </span><span class="kw">replicate</span>(<span class="dv">100</span>, <span class="kw">D</span>(<span class="kw">seq_len</span>(<span class="dv">128</span>)), <span class="dt">simplify =</span> <span class="ot">FALSE</span>)</span>
<span id="cb6-2"><a href="#cb6-2"></a>m_v &lt;-<span class="st"> </span><span class="kw">do.call</span>(cbind, <span class="kw">Map</span>(<span class="cf">function</span>(d) {d<span class="op">$</span>v}, x))</span>
<span id="cb6-3"><a href="#cb6-3"></a>m_w &lt;-<span class="st"> </span><span class="kw">do.call</span>(cbind, <span class="kw">Map</span>(<span class="cf">function</span>(d) {d<span class="op">$</span>w}, x))</span>
<span id="cb6-4"><a href="#cb6-4"></a></span>
<span id="cb6-5"><a href="#cb6-5"></a><span class="kw">library</span>(rhosa)</span>
<span id="cb6-6"><a href="#cb6-6"></a></span>
<span id="cb6-7"><a href="#cb6-7"></a>bc_v &lt;-<span class="st"> </span><span class="kw">bicoherence</span>(m_v, <span class="dt">window_function =</span> <span class="st">&#39;hamming&#39;</span>)</span>
<span id="cb6-8"><a href="#cb6-8"></a>bc_w &lt;-<span class="st"> </span><span class="kw">bicoherence</span>(m_w, <span class="dt">window_function =</span> <span class="st">&#39;hamming&#39;</span>)</span></code></pre></div>
<p>In the above code, we take 100 samples of the same length for a smoother result. The <code>bicoherence</code> function accepts a matrix whose column represents a sample sequence, and returns a data frame. Note that an optional argument to <code>bicoherence</code> is given for requesting tapering with Hamming <a href="https://en.wikipedia.org/wiki/Window_function">window function</a>.</p>
<div class="sourceCode" id="cb7"><pre class="sourceCode r"><code class="sourceCode r"><span id="cb7-1"><a href="#cb7-1"></a><span class="kw">library</span>(ggplot2)</span>
<span id="cb7-2"><a href="#cb7-2"></a></span>
<span id="cb7-3"><a href="#cb7-3"></a>plot_bicoherence &lt;-<span class="st"> </span><span class="cf">function</span>(bc) {</span>
<span id="cb7-4"><a href="#cb7-4"></a>    <span class="kw">ggplot</span>(bc, <span class="kw">aes</span>(f1, f2)) <span class="op">+</span></span>
<span id="cb7-5"><a href="#cb7-5"></a><span class="st">        </span><span class="kw">geom_raster</span>(<span class="kw">aes</span>(<span class="dt">fill =</span> value)) <span class="op">+</span></span>
<span id="cb7-6"><a href="#cb7-6"></a><span class="st">        </span><span class="kw">scale_fill_gradient</span>(<span class="dt">limits =</span> <span class="kw">c</span>(<span class="dv">0</span>, <span class="dv">10</span>)) <span class="op">+</span></span>
<span id="cb7-7"><a href="#cb7-7"></a><span class="st">        </span><span class="kw">coord_fixed</span>()</span>
<span id="cb7-8"><a href="#cb7-8"></a>}</span></code></pre></div>
<p>The axis <code>f1</code> and <code>f2</code> represent normalized frequencies in unit cycles/sample of range <code>[0, 1)</code>. Frequency pairs of bright points in the following plot of <code>bc_v</code> indicate the existence of some quadratic phase coupling, as expected:</p>
<div class="sourceCode" id="cb8"><pre class="sourceCode r"><code class="sourceCode r"><span id="cb8-1"><a href="#cb8-1"></a><span class="kw">plot_bicoherence</span>(bc_v)</span></code></pre></div>
<div class="figure">

<img src="" alt="v&#39;s estimated magnitude-squared bicoherence." width="50%" />

<p class="caption">

<p>v’s estimated magnitude-squared bicoherence.</p>
</p>

</div>

<p>In contrast, <code>bc_w</code> has no peaks at frequency pair <code>(f_1, f_2) = (0.35, 0.2)</code>, etc.:</p>
<div class="sourceCode" id="cb9"><pre class="sourceCode r"><code class="sourceCode r"><span id="cb9-1"><a href="#cb9-1"></a><span class="kw">plot_bicoherence</span>(bc_w)</span></code></pre></div>
<div class="figure">

<img src="" alt="w&#39;s estimated magnitude-squared bicoherence." width="50%" />

<p class="caption">

<p>w’s estimated magnitude-squared bicoherence.</p>
</p>

</div>

<h2 id="license">License</h2>
<p>GPLv3</p>
<h2 id="acknowledgement">Acknowledgement</h2>
<p>The author thanks Alessandro E. P. Villa for his generous support to this project.</p>
<h2 id="references">References</h2>
<p>[1] Brillinger, D. R., &amp; Irizarry, R. A. (1998). An investigation of the second- and higher-order spectra of music. Signal Processing, 65(2), 161–179. <a href="https://doi.org/10.1016/S0165-1684(97)00217-X">https://doi.org/10.1016/S0165-1684(97)00217-X</a></p>
<p>[2] Villa, A. E. P., &amp; Tetko, I. V. (2010). Cross-frequency coupling in mesiotemporal EEG recordings of epileptic patients. Journal of Physiology-Paris, 104(3), 197–202. <a href="https://doi.org/10.1016/j.jphysparis.2009.11.024">https://doi.org/10.1016/j.jphysparis.2009.11.024</a></p>

</body>
</html>
